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Three types of surface tensions can be defined for lipid membranes: the internal tension, a, conju¬ 
gated to the real membrane area in the Hamiltonian, the mechanical frame tension, r, conjugated to 
the projected area, and the “fluctuation tension”, r, obtained from the fluctuation spectrum of the 
membrane height. We investigate these surface tensions by means of a Monge gauge lattice Monte 
Carlo simulation involving the exact, nonlinear, Helfrich Hamiltonian and a measure correction for 
excess entropy of the Monge gauge. Our results for the relation between a and r agrees well with 
the theoretical prediction of [J.-B. Fournier and C. Barbetta, Phys. Rev. Lett., 2008, 100, 078103] 
based on a Gaussian approximation. This provides a valuable knowledge of r in the standard Gaus¬ 
sian models where the tension is controlled by a. However, contrary to the conjecture in the above 
paper, we find that r exhibits no significant difference from r over more than five decades of tension. 

Our results appear to be valid in the thermodynamic limit and are robust to changing the ensemble 
in which the membrane area is controlled. 


I. INTRODUCTION 

A few decades of research has clarified that lipid bilay¬ 
ers can be treated as sheets obeying the continuum curva¬ 
ture/tension elasticity formulated by Helfrich]!]. Char¬ 
acterising membrane elastic constants in numerical sim¬ 
ulations has become an important issue due to the in¬ 
crease in computational power and development of sim¬ 
ulation methods. Starting with the seminal work of 
Ref. 2, a growing number of research papers are dedicated 
to the determination of bending rigidity [3-12], Gaus¬ 
sian modulus[13], and surface tension]?, 8, 11, 14, 15], 
from molecular simulation of all-atom systems to coarse¬ 
grained models. 

A number of questions have been raised regarding the 
definition of membrane surface tension. When a ther¬ 
mally fluctuating flexible membrane is stretched by a 
moderate lateral tension, the membrane extends its pro¬ 
jected area by suppressing its transverse fluctuation[16]. 
This leads to multiple definitions of membrane tension: 
i) “internal tension”, a, conjugated to the real area A in 
the membrane Hamiltonian, ii) mechanical “frame ten¬ 
sion”, T, conjugated to the projected area Ap, [17] and 
iii) “fluctuation tension”, r, associated to the lowest- 
order wavevector dependence of the inverse fluctuation 
spectrum. The latter two are directly experimentally 
measurable[18]. 

Some discussions on the difference between these sur¬ 
face tensions have been lasting for more than a decade. 
First of all, it is well accepted that a and t are in¬ 
trinsically different, because cr is a microscopic param¬ 
eter belonging to the effective Hamiltonian while r is a 
macroscopic observable integrating the local stress ten¬ 
sor over all the membrane fluctuations. Based on free- 
energy calculations and stress tensor averages, the dif¬ 


ference between a and r has been estimated within the 
Gaussian approximation [14, 19, 20]. Note that the va¬ 
lidity of this estimation was questioned by Schmid [21]. 
It is also well accepted that a and r are different, be¬ 
cause r is actually the renormalised version of a (at the 
Gaussian level they are equal). The relationship between 
r and r has been much debated. Several authors have 
argued theoretically[5, 15, 22, 23], or observed in numer¬ 
ical simulations[2-4, 7, 8, 15, 21, 24-28], that r matches 
the mechanical frame tension r. Other authors have ob¬ 
served a difference between r and r, either in numerical 
simulations [11, 14, 20, 29] or experiments[18]. 

These questions should be carefully revisited with the 
help of extensive, large-scale numerical simulations. A 
complication has to do with the equivalence of ensem¬ 
bles in the thermodynamic sense, as in the simulations 
found in the literature, either the areas A and Ap or their 
conjugated tensions a and r are fixed. Different ensem¬ 
bles are only equivalent in the thermodynamic limit of 
infinite membranes, A —>■ oo at hxed A/N, with N the 
total number of degrees of freedom. It is therefore very 
important to check the convergence of the numerical re¬ 
sults with system size. In another recent simulation[30], 
r was claimed to deviate from r for t < 0 and to coin¬ 
cide with r for T > 0. It is thus interesting to study how 
r behaves when r assumes positive values very close to 
zero, for which the system undergoes large fluctuations. 
Let us stress that investigating several orders of magni¬ 
tude of r is important, because membrane tension can 
experimentally vary over more than five decades[31]. 

In this paper, we investigate the differences between 
the three surface tensions, a, t and r, over several 
decades of t, by means of Monte Garlo simulations of a 
lattice membrane system in the Monge gauge. In the lat¬ 
tice model, the internal tension cr and the frame tension 
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T can be controlled and short range fluctuations gener¬ 
ated by molecular or particle protrusions, which make 
accurate estimations of r difficult, are not present. How¬ 
ever, a measure correction is required to cure the excess 
entropy due to the unidirectional lattice-site motions im¬ 
plied by the Monge gauge. We also employ exact non¬ 
linear expressions to calculate the membrane’s curvature 
and area. We find that cr differs from r at small positive 
frame tensions, while they match at moderate and high 
frame tensions. Our results for a — t agree well with the 
prediction of ref. 20 based on a Gaussian approximation, 
although our simulation investigates regimes where the 
Gaussian approximation is far from being justified. We 
find also that r and t match, within the accuracy of our 
simulation, over more than five decades, even in the limit 
of very small frame tensions. These relationships remain 
in the thermodynamic limit and seem not to depend on 
the ensemble in which the real area A is controlled (see 
Appendix IV A). 

Our paper is organised as follows. In Sec. II, we detail 
our lattice membrane model and the Monte Garlo sim¬ 
ulation framework that realises an equilibrium ensemble 
for fixed r and a. In this framework, the average of the 
real membrane area, (A), is kept constant by adjusting cr 
while T is varied. We also perform our simulation in the 
equilibrium ensemble for fixed r and A. This is discussed 
in Appendix IV A. In the thermodynamic limit of large 
systems, we thus simulate a membrane with fixed real 
area subjected to a variable frame tension. In Sec. Ill, 
we numerically determine and we compare the three ten¬ 
sions cr, r and r. We give our conclusions in Sec. IV. 


II. MODEL AND METHOD 

A. Stress-controlled ensemble 

Let us consider an incompressible, fluctuating lipid 
membrane with a fixed area A = Aq, corresponding to a 
fixed number of lipids. We assume that it is attached to 
a deformable frame, of variable area Ap, that exerts onto 
the membrane a fixed tension t. The Hamiltonian of the 
system is thus given by[I] 

Ti =—rAp + I dA, (r, A fixed; Ap free) (1) 

Ja 2 

where H = ci + C 2 represents the sum of the two prin¬ 
cipal curvatures of the membrane. The Gaussian bend¬ 
ing modulus [I] R has been neglected, according to the 
Gauss-Bonnet theorem, which is correct for fixed angu¬ 
lar boundary conditions or periodic boundary conditions 
along the frame [32]. 

We use the Monge representation, in which the mem¬ 
brane shape is describe by its height z = h{x,y) above 
the plane of the frame. Note that this parametrisation 
does not require the membrane deformation to be small. 
The only restriction is that overhangs are forbidden, since 


the function h is single valued. In this representation, the 
total curvature H is exactly given by[33] 

(1 -|- h‘^)hyy -|- (1 -f hy)hxx ~ ‘^hxhyhxy 

- il + hl + hlfA ’ ^ ^ 

where subscripts represent spatial derivatives, e.g., hx = 
dhjdx. We work here in the non-linear regime, allow¬ 
ing for large deformations of the membrane, i.e., we do 
not use the custom linear approximation H ~ hxx + hyy 
(employed in the Gaussian framework). Note that in our 
simulations (see the details below) the maximum slope, 
max(|Ii 2 ;|), reaches 1.3 at vanishing frame tension, but 
does not exceed this value. This justifies both using the 
exact expression for H and neglecting overhangs. Close 
to the buckling transition, however, i.e., for large nega¬ 
tive frame tensions, these assumptions may not be valid. 

Because it is quite difficult to work with a fixed area 
A = Aq, we change ensemble in order to control the 
conjugated internal tension cr, instead of the membrane 
area A. In the thermodynamic limit, the two ensembles 
are expected to be equivalent. Indeed, as shown in Ap¬ 
pendix IV A, we find quantitatively similar results in the 
ensemble in which the real area A is almost prescribed by 
means of a quadratic potential. In this stress-controlled 
ensemble, the modified Hamiltonian weighting the mem¬ 
brane fluctuations is given by 

'H* = —tAp+ ^11“^ dA+aA, (r,cr fixed; Ap,Afree) 
Ja 2 

( 3 ) 

with (A) = Aq systematically enforced by adjusting a. 

B. Lattice model 

In order to implement the membrane fluctuations nu¬ 
merically, we introduce a Nx x Ny lattice on the plane 
{x, y) of the frame, with lattice spacing a. The height of 
the membrane surface is defined by hij on each lattice 
site, with 1 < z < Nx and 1 < j < Ny. In the reminder 
of this paper, we take Nx = Ny = '/N. Furthermore, we 
adopt periodic boundary condition. Thus the projected 
area is given by Ap = Va^. 

Importantly, because the number of lipids is fixed while 
the frame area Ap is variable, we take N constant but a 
variable. In other words, in our simulation the lattice 
spacing is not fixed. We stress, however, that a always 
remains uniform over the lattice: it is the overall lattice 
spacing that changes. Therefore the projection over the 
reference plane of the vertices always maps onto a square 
grid of finite spacing and the vertices can never overlap. 
Note that a similar lattice membrane model (with fixed 
spacing and the Gaussian approximation) has been em¬ 
ployed for Monte Garlo simulation in order to investigate 
inclusions effects [34]. 

Monte Garlo simulations are performed with the lat¬ 
tice heights {hij} and the lattice spacing a as dynami¬ 
cal variables. Special care is required, however, for the 
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sampling of the height variables. In the Monge gauge, 
each site can only fluctuate in the z direction, whereas in 
the physical gauge, because the membrane is fluid, the 
relevant fluctuations occur by local displacements along 
the membrane’s normal. In order to avoid an artificial 
entropy increase of the states where the membrane is 
tilted, the naive measure d^[h] = ^ dhij/So, where Jq 

is a quantum of height fluctuations, must be transformed 
to 


n 


dhij 

do/ cos 9ij ■ 


(4) 



80 


Here, 9ij is the angle between the normal vector of the 
site (f,j) and the z axis, and cos 9ij = 1/.^! + + /i^. 

In this way, the quantum height displacement in the nor¬ 
mal direction is always Sq. Equivalently, we can keep 
the naive measure d^ [h] and add to the Hamiltonian the 
following correction term: 

■Hcorr = -knT ^ln(cos6»^j), (5) 


FIG. 1. Simulation snapshot for k = 10 ksT, r = 0 (ten¬ 
sionless frame), N = 6400, and {A) — Nuq achieved for 
a ~ 0.50 ksT/uQ. The lengths are given in units of no, 
which corresponds to the size of a coarse-grained degree of 
freedom in the membrane. In this simulation, the equilib¬ 
rium lattice spacing is a ~ 0.983 ao, which corresponds to 
{A — Ap)/A ~ 0.034. With typically ao ~ 20 nm and 
/cbT = 4 X 10~^^ J, the residual internal tension is about 
cr ~ 5 X 10“® J/m'^. 


where T is the temperature and kB Boltzmann’s con¬ 
stant. Note that this correction, that would be exact 
if the membrane was uniformly tilted, is only approxi¬ 
mate if the tilt changes from site to site. In practice, it 
works well when the membrane slope is small, as in the 
present simulations. For highly tilted membranes, such 
as buckled membranes, the correction term, eqn (5), is 
not sufficient, as discussed in Appendix IV B. 

The partition function in the stress-controlled ensem¬ 
ble introduced above is given by 

ZN{T,a) = J J J d^[h] exp[-H7(fcBT)], (6) 

with %' = %* + Hcorr- Metropolis sampling is em¬ 
ployed [35]. In each Monte Carlo step, N trial moves of 
the height of a random site are attempted. The change 
of a (hence Ap) is attempted every 5 Monte Carlo steps. 
The amplitudes of the hij and Ap changes are adjusted so 
that the rejection rate lies in the range 40%-60%. Both 
A and Ap change during the simulations. However, for a 
given Monte Carlo run at a specified value of r, the inter¬ 
nal tension cr is adjusted in order for the average of the 
membrane area (A) to have the specihc value Ag, fixed 
once for all. 

In our simulation, the real area A of the membrane is 
calculated as follows. The local area associated to a site 
P{i,j) is obtained as one half of the total area of the four 
triangles [APB], [BPC], [CPD] and [DPA] built from the 
neighbouring sites of coordinates: 

A(j + 1, j). B(z,j -b 1), C(f - 1, j), D(f, J - !)• (7) 

The area of each triangle can be calculated by using 
height differences between the lattice sites. Then, the 


total real area can be expressed as 

NN 

^ = EE x(Tapb + Tbpc + Tcpd + Tdpa)- (8) 

i=i j=i ^ 

Simulation are carried out at least for 10^ steps after 
the initial relaxation. The statistical errors are calculated 
from four or more independent runs. 

C. Dimensionless units and orders of magnitudes 

Since our simulated membrane has a prescribed aver¬ 
age area Ag and a hxed number V of degrees of freedom, 
we can picture each degree of freedom as a patch of lipids 
of fixed average area a§ = Ag/V. In the following we are 
going to nondimensionalise the lengths by ag (the size 
of a degree of freedom in the membrane) and the ener¬ 
gies by /cpT = 4 X 10“^^ J (room temperature energy). 
In other words, we shall set og = /cbT = 1. Hence, all 
tensions will be given in units of feeT/a^ As previously 
discussed, we shall thus adjust a for each r in order to 
have systematically (A) = V, in dimensionless units. 

Let us estimate tension scale feeT/ag. In principle, we 
are free to choose the linear size ag of the lipid patches 
that constitute our degrees of freedom. It is convenient, 
however, to choose the smallest possible size, in order 
to allow for all possible fluctuations and to avoid us¬ 
ing renormalised elastic constants. Owing to the level 
of coarse-graining of our simulation, in which the bi¬ 
layer membrane of real thickness ~ 5 nm is treated as 
a mathematical surface, it is natural to think of the lipid 
patches as regions of typical size, e.g. ag ~ 20 nm. 
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Then, the natural unit of tension in our simulation is 
k-gT/a^ K, 10“® J/m^. The tension mechanically im¬ 
posed on ordinary vesicle membranes, in experiments, are 
in the range 10“^ J/m^-10“^ J/ni^, where the smallest 
value corresponds to floppy vesicles and the largest value 
corresponds to the lytic tension[36]. We shall therefore 
have r span the range 10“^-10^ in dimensionless units. 
Let us stress that this corresponds to five orders of mag¬ 
nitude. 

In our simulation, we choose a specific value of the 
bending rigidity corresponding to k = 10 /cbT, i.e., k = 
10 in dimensionless units. We shall consider membrane 
sizes corresponding to = 400,1600,6400 and 25600. 
With the value of ag given above, this corresponds to 
membrane of linear dimension \/Aq = ^/Nag = 400 nm, 
800 nm, 1.6 /im and 3.2 /rm, respectively. These values 
are quite small if one thinks of macroscopic experiments, 
but they are reasonably large on the biological scale. Be¬ 
cause of the change of ensemble that we have performed, 
we shall check the convergence of our results with in¬ 
creasing N. 

In Fig. 1, a simulation snapshot for N = 6400 is shown 
in the tensionless state r = 0. 


III. RESULTS 

A. Frame tension and internal tension 



FIG. 2. (a) The values of a corresponding to the prescribed 
average membrane area (A) = N are plotted against |t|, for 
system sizes N = 400,1600, 6400 and 25600. The upper (resp. 
lower) branch displays the data for r > 0 (resp. r < 0). 
For r < 0 the data points are plotted up to the buckling 
transition. The dotted line indicates a = t (for r > 0). (b) 
The red points show the difference ct —r as a function of r > 0 
for N = 6400; the black dashed line shows the best fit to the 
theoretical Gaussian prediction eqn (9)-(10). 


From now, unless otherwise specified, all quantities will 
be given in dimensionless units (see Sec. IIC). We investi¬ 
gate here the relationship between tr and r. As explained 
above, we place ourselves in the (cr, r) ensemble; however, 
for every imposed frame tension r, we determine the in¬ 
ternal tension cr that achieves (A) = N. Our simulated 
membrane thus effectively has a constant real area (up 
to fluctuations that become irrelevant in the thermody¬ 
namic limit). Simulations in the (A,r) ensemble, where 
A is nearly fixed by means of a steep quadratic potential, 
are discussed in Appendix IV A. 

The determined values of cr are plotted in Fig. 2a 
against the frame tension for different system sizes: N = 
400,1600, 6400 and 25600. For any value of r we find 
a > T. Data for both positive and negative frame ten¬ 
sions are displayed together as a function of |r|: the up¬ 
per branch corresponds to r > 0 and the lower branch 
corresponds to r < 0. Note that for r > 0 the difference 
between cr and r is sizeable only for small values of the 
frame tension. For r = 0 we find that the residual inter¬ 
nal tension is equal to ag ~ 0.500±0.001, where the error 
bar takes into account the determinations using different 
values of N. In other words, the residual tension turns 
out to be fairly independent of system size. 

When T decreases below some negative threshold < 
0, the membrane buckles either into the x or the y direc¬ 
tion: the average membrane shape undergoes a transition 
from flat to non-flat through a symmetry breaking (like a 
rod under compression). Note that the membrane height 


function remains single valued although it acquires a bi- 
modal distribution. The absolute value \ti,\ of the buck¬ 
ling threshold decreases as the system size N gets larger. 
In the lower branch, the data points for |t| > |Tt,|, in the 
buckling state, are eliminated from the plots in Figs. 2a 
and 3. It is notable that the condition cr < 0 is not re¬ 
quired for the buckling transition. We also simulated our 
lattice model in the (Ap,CT) ensemble, i.e. with constant 
lattice spacing a and freely changing membrane area A. 
We found that the membrane buckles then at tr = 0. 
These results suggest that cr is not the mechanical force 
that drives membrane buckling. 

In Ref.20, the difference between the tensions cr and r 
was calculated within the Gaussian approximation, yield¬ 
ing 


cr — T = 


k^Th^ 
Stt 


cr / /tA^ 

I- — In 1 -I- 

kA^ V cr 


( 9 ) 


where A is the ultraviolet wavevector cutoff and all quan¬ 
tities have their normal dimensions. Within our simu¬ 
lation, this ultraviolet cutoff corresponds to the lattice 
spacing. We thus expect 

A = Q! — , with cr ~ 1. (10) 

Whereas our simulation is fully nonlinear, the expression 
for cr — r calculated within the Gaussian approximation 
fits the data quite well (Fig. 2b). Since the numerical 
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(A-Ap)/A 


FIG. 3. The internal tension a (closed circles) and the frame 
tension r (open circles) are plotted as a function of residual 
area {A — Ap) j A for various system sizes N = 400,1600, 6400 
and 25600. The data for r is not shown in the region r < 0. 
Finite size effects becomes small as the system size increases: 
the data suggest that the thermodynamic limit is almost 
reached for N = 6400. All the statistical errors are small, 
within the symbol marks. 


value of (a) changes by only 1.4% from t = 0 to 1000, we 
use for the fit the extreme value at high tension, a = 1. 
The best fit of the data gives then a = 1.12, very close 
to unity. 

Thermodynamic limit .—In Fig. 3, we plot the two ten¬ 
sions cr and T against the excess area {A — Ap)/A for 
the various system sizes N = 400,1600, 6400 and 25600. 
Here, only positive values of t are displayed: the down¬ 
wards divergence of the curves associated to the plain 
circles corresponds to r —>■ 0. The system size depen¬ 
dence is marked for small sizes, but the data appears to 
converge to a universal thermodynamic limit at large sys¬ 
tem sizes. Finite size effects are more or less irrelevant 
for N > 6400. 


B. Fluctuation tension 


We now discuss the fluctuation tension r. It is derived 
from the spectrum density, assumed to take the standard 
form [33]: 


(IM9)P) 


-I- ’ 


( 11 ) 


where h{q) is the Fourier transform of the height func¬ 
tion, with q the wavevector. The fluctuation tension r 
and the bending rigidity Kr, appearing in the fluctua¬ 
tion spectrum, are the renormalised counterparts of the 
Hamiltonian parameters a and n. 



(A-Ap)/A 



FIG. 4. (a) Plot of r, r, and cr as a function of excess area 
{A — Ap)/A for N = 6400, which roughly corresponds to the 
thermodynamic limit. The data for r and cr are the same as in 
Fig. 3, and r is plotted with the error bars. As for the statis¬ 
tical errors in {A — Ap)/A, they are small, within the symbol 
marks, (b) The difference between r and r, normalised by 
CT, is plotted against r for N = 1600, 6400, and 25600. The 
convergence of the data with increasing N indicates that the 
thermodynamic limit is reached for N — 6400. For the sake 
of comparison, a plot of {a — T)/a is also shown. 


In estimating these quantities using simulation results, 
it is crucially important to consider large-scale fluctu¬ 
ations. However, because we always treat finite-size 
systems, the wavevector spectrum is limited, which re¬ 
sults in systematic deviations. Previously, two of the 
present authors have been addressed in detail the ques¬ 
tion of the estimation for the bending rigidity from simu¬ 
lation data[8]. For a planar membrane, the correct value 
of the bending rigidity can be determined uniquely by 
first estimating the rigidity Kr via a fit of the inverse 
spectrum {rq^ + Krq"^)/{k-eT) in the wavevector range 
2'K/{ao'/N) < q < (7cut, then by extrapolating the re¬ 
sults to the limit q^t —t 0, where ^cut is a varying upper 
fitting wavevector limit. Note that because the projected 
area Ap fluctuates, there is some ambiguity regarding the 
definition of the quantified wavevectors. As discussed in 
Appendix IV C, our results are not sensitive to those de¬ 
tails. This procedure has been carried out for the present 
system. We find that the value of Kr coincides with the 
bare value k, up to our numerical precision. This result 
makes sense because renormalisation group calculations 
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predict 

k T 

Ak = Kr — K -^( 12 ) 

Stt ’ ^ ^ 

up to a numerical factor of order unity[33] (here all quan¬ 
tities have their normal dimensions). This prediction 
gives Ak/k ~ 0.03 for N = 6400, which is less than 
the error bars, and thus not detectable in our simulation. 

The correct value of the fluctuation tension r (the 
renormalised tension) can be estimated likewise in the 
limit gcut —0. For finite values of ^cut, we find that r 
is overestimated but converges in the limit q^nt —>■ 0 to 
a well defined value. The corresponding extrapolation 
is performed by using a quadratic function of (/cut- The 
fitting was performed in the range 0.015 < [qout/'^Y < 
0.15. 

Figure 4 shows the values of the fluctuation tension r. 
We find that r is fairly close to r for all the investigated 
values of r, as indicated in Fig. 4a. In Fig. 4b, the 
difference r — r is plotted against r. For large system 
sizes, N = 6400 and 25600, we observe a small deviation 
between r and r; however, as discussed in Appendix IV B, 
we believe that this deviation is smaller than a possible 
systematic error arising from the measure correction. 


C. Comparison with the Gaussian approximation 

In the present work, all the nonlinearities of the prob¬ 
lem have been taken into account: (i) the nonlinear ex¬ 
pression of the membrane elementary area dA (ii) the 
nonlinear expression of the membrane bending energy 
density H^dA, and (iii) the nonlinear measure correction 
"Hcorr dealing with the excess entropy of the Monge gauge. 
The Gaussian approximation, which is frequently used in 
basic calculations, consists in replacing the bending en¬ 
ergy density H'^dA by HldAp, where Hi, = h^x + hyy is 
the Laplacian approximation of the mean-curvature and 
dAp the projected elementary area, and in neglecting the 
entropy measure correction "Hcorr- Indeed, at quadratic 
order, the measure correction becomes ~ 

which is of no consequences as it sim¬ 
ply redefines the tunable parameter cr. In this subsection, 
we investigate how these two Gaussian approximations 
affect the surface tensions. First, we make the Lapla¬ 
cian approximation, i.e., H^dA is replaced by H^dAp in 
eqn (1), and we keep the entropy correction term Hcoir- 
The corresponding results are denoted by the subscript 
EL in Fig. 5. Then we also remove the entropy correction 
term. The corresponding calculations are denoted by the 
subscript L. 

We find that the internal tension cr is not changed by 
these modifications. However, tel and rn are signifi¬ 
cantly different from r, as shown in Fig. 5. Using the 
Laplacian approximation while keeping the entropy cor¬ 
rection term results in a fluctuation tension tel about 
three times larger than a. Further removing the entropy 



FIG. 5. Excess fluctuation tensions r, tel and tl, measured 
with respect to r and normalized by cr, plotted against r, 
for the nonlinear case (r), and for the Gaussian curvature 
approximation with the entropy correction (rBc) or without 
it (ru). Here, N = 6400. For comparison, a plot of (cr — r)/cr 
is shown. 


correction yields a tension tl very close to cr, in agree¬ 
ment with the well-known coincidence of r and a at the 
Gaussian level. Therefore, the effects of the measure cor¬ 
rection and of the nonlinear curvature are opposite to 
each other: the former raises the value of r, the latter 
lowers it. Both of these effects need to be taken into 
account for a faithful description of the membrane fluc¬ 
tuations beyond the level of the Gaussian approximation, 
and the estimated r deviates from the true value if either 
of them is missing. Note that without the entropy cor¬ 
rection term, Hcorr, the nonlinear bending energy leads 
to an enhancement of the vertical motion of the sites, 
resulting in the collapse of the simulations. 


IV. CONCLUSIONS 

In this paper, we have examined the various surface 
tensions of membranes: frame tension r, internal ten¬ 
sion cr, and renormalised “fluctuation” tension r. We 
have compared them quantitatively using a Monte-Carlo 
simulation with control over both r and cr, the latter 
being slaved to the former in order to keep the average 
membrane area (A) constant. We have also validated 
our results in the conjugated ensemble where the real 
membrane area A is fixed (Appendix IV A). We have in¬ 
vestigated large systems in order to reach the thermody¬ 
namic limit. Our model being a lattice model instead of 
a particle-based model, protrusion effects and other ar- 
tihcial elements are eliminated. Gaussian approximation 
artefacts are also excluded because curvatures and areas 
are computed with exact formulas (although we employ 
the Monge gauge) and because we have corrected the ex¬ 
cess entropy associated with the measure of the Monge 
gauge. Moreover, because our simulated membrane does 
not exhibit rupture at extremely large r, all tensions from 
vanishingly small to very high can be investigated. 
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There have been two findings in our results. The first 
is that the residual internal tension ao = a — t remains 
finite at large N (large systems). While it is well known 
that cr and r are different, we have confirmed that the 
theoretical (albeit Gaussian) estimation of Ref.20 pre¬ 
dicts correctly and quantitatively the difference a — t 
from vanishingly small frame tensions to very large frame 
tensions. 

The renormalised “fluctuation” tension r has also been 
investigated, and compared with the other two tensions. 
We conclude r ~ t within an accuracy of Q.lk-QT/a^ 
(estimated ~ 10“® J/m ), which is consistent with the 
proposition of Ref.22. 

Note finally that the shape of the buckled membrane, 
described by elliptic functions, breaks the symmetry be¬ 
tween the X and y directions and yields an anisotropy in 
the frame tension [9]. How exactly the thermal fluctua¬ 
tions of the membrane modify the mechanical buckling 
transition should be investigated in further studies. 


APPENDIX 

A. Simulation with nearly constant real area 

In the body of this paper, the real area A of the mem¬ 
brane is not fixed, but its average value (A) is controlled 
by the parameter a (much as the average number of par¬ 
ticles in a grand-canonical ensemble is controlled by the 
chemical potential). In principle, we have checked the 
validity our results by studying the thermodynamic limit 
of increasing system sizes. It is nonetheless interesting to 
test our results by working in the ensemble in which the 
area A is prescribed. This is the purpose of the present 
Appendix. 

While it is very difficult to prescribe exactly the mem¬ 
brane area, because it is difficult to identify the sampling 
condition that satisfy this constraint, one can easily al¬ 
most prescribe the area by using a quadratic potential of 
tunable strength Ka- We thus worked with the Hamil¬ 
tonian 

n* = -rAp + [ dA + ^{A - Ao)2. (13) 



FIG. 6. The deviation of renormalised tension r from the 
frame tension r is plotted in the form (r — x^jo. The red 
points are pasted from Fig. 4b for comparison. The results 
for three values of Ka (=1, 10, and 100) are shown by using 
symbols and error bars of different colours. For all the plotted 
data, the system size is A = 6400. Although we normalise the 
data by cr, no estimate is available in the present simulation; 
therefore, we use the values of Fig. 4a as substitutes. 


B. Entropy correction for membrane tilt 

To check the reliability of the measure correction "Hcorr, 
eqn (5), we performed simulations in which the en¬ 
tire membrane is tilted by an angle 0^. The periodic 
boundary condition along the x direction is modified to 
h{Nx,j) = h{0,j) + N^asinOa- The x side length of 
the lattice is changed to a cos 9a in order to maintain a 
square grid in the tilted projected plane. The height spec¬ 
trum is calculated for the tilted plane, of normal vector 
(—sin 0a, 0, cos da)- The fluctuation tension r increases 
with increasing 9a (see Fig. 7), while tr is independent of 
da- This implies that the excess entropy associated with 
the membrane tilt in the Monge gauge is not completely 
removed. However, the deviation is small for small da, 
and the overestimation of r is less than 0.1 for da < O.Itt. 
In our simulation, at r = 0 we obtain (cosdy) = 0.98, 
i.e., the mean angle is O.OOOtt. We conclude that the ob¬ 
tained values of r are reliable within the accuracy of 0.1 
(in dimensionless units). 


C. Definition of the spectrum wavevectors 


We have performed the corresponding Monte-Carlo 
simulations exactly in the same manner as in the body of 
the paper. Figure 6 shows the result of these additional 
simulations, for Ka = 1, 10 and 100, with N = 6400, and 
compares them with the results of Fig. 4. With increas¬ 
ing KA, the standard deviation of the real area decreases: 
0.01%, 0.001%, and 0.0003%, for Ka = 1, 10, and 100, 
respectively. With larger values of Ka, the fluctuations 
of the membrane area become suppressed so that unbi¬ 
ased equilibrium sampling is impossible. 


In Sec. HIB, the renormalised tension r is calculated 
by fitting the height spectrum (|/i(q)p). We use a se¬ 
quence of instantaneous height configurations, and we 
perform averages after the simulation finishes. Because 
the projected area Ap changes during the simulations, 
the wavevectors on the reciprocal lattice are not fixed. 
They are given by 


9 — (da;, dy) — 


27r 

\/^ 


(n,m), 


(14) 















FIG. 7. Dependence of the fluctuation tension r for a mem¬ 
brane whose projected plane is entirely tilted by the angle 0a, 
for T — Q and r = 1, at = 6400. 

where n and m are integers. All these data are mapped 
to a one-dimensional wavevector q = \q\. In plotting 
(|/i(g)p), the data of \h{q)\^ are averaged over close val¬ 
ues of q (the width of q bins is bir/'/N. Here, 6 = 1.2 is 
taken) and over all the height configurations. 

We have confirmed the robustness of our estimation of 
r, by using a different set of wavevectors, defined by 

/ 27F . , , . 

q=^=={n,m). (15 

vW) 


For each of integers n,m, the values of {\h{q')\'^) is aver¬ 
aged over the sequence of instantaneous height configu¬ 
rations. The results obtained by these two methods do 
not essentially differ from one another: the differences in 
the value of r are smaller than 0.025 k^T/aQ for all r, 
which is typically smaller than the error bars. 
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